00001
00002
00003
00004 #include <stdio.h>
00005 #include "graph.h"
00006 #include "instances.inc"
00007
00008
00009
00010
00011
00012 #define TERMINAL ( (arc *) 1 )
00013 #define ORPHAN ( (arc *) 2 )
00014
00015
00016 #define INFINITE_D ((int)(((unsigned)-1)/2))
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034 template <typename captype, typename tcaptype, typename flowtype>
00035 inline void Graph<captype,tcaptype,flowtype>::set_active(node *i)
00036 {
00037 if (!i->next)
00038 {
00039
00040 if (queue_last[1]) queue_last[1] -> next = i;
00041 else queue_first[1] = i;
00042 queue_last[1] = i;
00043 i -> next = i;
00044 }
00045 }
00046
00047
00048
00049
00050
00051
00052 template <typename captype, typename tcaptype, typename flowtype>
00053 inline typename Graph<captype,tcaptype,flowtype>::node* Graph<captype,tcaptype,flowtype>::next_active()
00054 {
00055 node *i;
00056
00057 while ( 1 )
00058 {
00059 if (!(i=queue_first[0]))
00060 {
00061 queue_first[0] = i = queue_first[1];
00062 queue_last[0] = queue_last[1];
00063 queue_first[1] = NULL;
00064 queue_last[1] = NULL;
00065 if (!i) return NULL;
00066 }
00067
00068
00069 if (i->next == i) queue_first[0] = queue_last[0] = NULL;
00070 else queue_first[0] = i -> next;
00071 i -> next = NULL;
00072
00073
00074 if (i->parent) return i;
00075 }
00076 }
00077
00078
00079
00080 template <typename captype, typename tcaptype, typename flowtype>
00081 inline void Graph<captype,tcaptype,flowtype>::set_orphan_front(node *i)
00082 {
00083 nodeptr *np;
00084 i -> parent = ORPHAN;
00085 np = nodeptr_block -> New();
00086 np -> ptr = i;
00087 np -> next = orphan_first;
00088 orphan_first = np;
00089 }
00090
00091 template <typename captype, typename tcaptype, typename flowtype>
00092 inline void Graph<captype,tcaptype,flowtype>::set_orphan_rear(node *i)
00093 {
00094 nodeptr *np;
00095 i -> parent = ORPHAN;
00096 np = nodeptr_block -> New();
00097 np -> ptr = i;
00098 if (orphan_last) orphan_last -> next = np;
00099 else orphan_first = np;
00100 orphan_last = np;
00101 np -> next = NULL;
00102 }
00103
00104
00105
00106 template <typename captype, typename tcaptype, typename flowtype>
00107 inline void Graph<captype,tcaptype,flowtype>::add_to_changed_list(node *i)
00108 {
00109 if (changed_list && !i->is_in_changed_list)
00110 {
00111 node_id* ptr = changed_list->New();
00112 *ptr = (node_id)(i - nodes);
00113 i->is_in_changed_list = true;
00114 }
00115 }
00116
00117
00118
00119 template <typename captype, typename tcaptype, typename flowtype>
00120 void Graph<captype,tcaptype,flowtype>::maxflow_init()
00121 {
00122 node *i;
00123
00124 queue_first[0] = queue_last[0] = NULL;
00125 queue_first[1] = queue_last[1] = NULL;
00126 orphan_first = NULL;
00127
00128 TIME = 0;
00129
00130 for (i=nodes; i<node_last; i++)
00131 {
00132 i -> next = NULL;
00133 i -> is_marked = 0;
00134 i -> is_in_changed_list = 0;
00135 i -> TS = TIME;
00136 if (i->tr_cap > 0)
00137 {
00138
00139 i -> is_sink = 0;
00140 i -> parent = TERMINAL;
00141 set_active(i);
00142 i -> DIST = 1;
00143 }
00144 else if (i->tr_cap < 0)
00145 {
00146
00147 i -> is_sink = 1;
00148 i -> parent = TERMINAL;
00149 set_active(i);
00150 i -> DIST = 1;
00151 }
00152 else
00153 {
00154 i -> parent = NULL;
00155 }
00156 }
00157 }
00158
00159 template <typename captype, typename tcaptype, typename flowtype>
00160 void Graph<captype,tcaptype,flowtype>::maxflow_reuse_trees_init()
00161 {
00162 node* i;
00163 node* j;
00164 node* queue = queue_first[1];
00165 arc* a;
00166 nodeptr* np;
00167
00168 queue_first[0] = queue_last[0] = NULL;
00169 queue_first[1] = queue_last[1] = NULL;
00170 orphan_first = orphan_last = NULL;
00171
00172 TIME ++;
00173
00174 while ((i=queue))
00175 {
00176 queue = i->next;
00177 if (queue == i) queue = NULL;
00178 i->next = NULL;
00179 i->is_marked = 0;
00180 set_active(i);
00181
00182 if (i->tr_cap == 0)
00183 {
00184 if (i->parent) set_orphan_rear(i);
00185 continue;
00186 }
00187
00188 if (i->tr_cap > 0)
00189 {
00190 if (!i->parent || i->is_sink)
00191 {
00192 i->is_sink = 0;
00193 for (a=i->first; a; a=a->next)
00194 {
00195 j = a->head;
00196 if (!j->is_marked)
00197 {
00198 if (j->parent == a->sister) set_orphan_rear(j);
00199 if (j->parent && j->is_sink && a->r_cap > 0) set_active(j);
00200 }
00201 }
00202 add_to_changed_list(i);
00203 }
00204 }
00205 else
00206 {
00207 if (!i->parent || !i->is_sink)
00208 {
00209 i->is_sink = 1;
00210 for (a=i->first; a; a=a->next)
00211 {
00212 j = a->head;
00213 if (!j->is_marked)
00214 {
00215 if (j->parent == a->sister) set_orphan_rear(j);
00216 if (j->parent && !j->is_sink && a->sister->r_cap > 0) set_active(j);
00217 }
00218 }
00219 add_to_changed_list(i);
00220 }
00221 }
00222 i->parent = TERMINAL;
00223 i -> TS = TIME;
00224 i -> DIST = 1;
00225 }
00226
00227
00228
00229
00230 while ((np=orphan_first))
00231 {
00232 orphan_first = np -> next;
00233 i = np -> ptr;
00234 nodeptr_block -> Delete(np);
00235 if (!orphan_first) orphan_last = NULL;
00236 if (i->is_sink) process_sink_orphan(i);
00237 else process_source_orphan(i);
00238 }
00239
00240
00241
00242 }
00243
00244 template <typename captype, typename tcaptype, typename flowtype>
00245 void Graph<captype,tcaptype,flowtype>::augment(arc *middle_arc)
00246 {
00247 node *i;
00248 arc *a;
00249 tcaptype bottleneck;
00250
00251
00252
00253
00254 bottleneck = middle_arc -> r_cap;
00255 for (i=middle_arc->sister->head; ; i=a->head)
00256 {
00257 a = i -> parent;
00258 if (a == TERMINAL) break;
00259 if (bottleneck > a->sister->r_cap) bottleneck = a -> sister -> r_cap;
00260 }
00261 if (bottleneck > i->tr_cap) bottleneck = i -> tr_cap;
00262
00263 for (i=middle_arc->head; ; i=a->head)
00264 {
00265 a = i -> parent;
00266 if (a == TERMINAL) break;
00267 if (bottleneck > a->r_cap) bottleneck = a -> r_cap;
00268 }
00269 if (bottleneck > - i->tr_cap) bottleneck = - i -> tr_cap;
00270
00271
00272
00273
00274 save(middle_arc -> sister);
00275 middle_arc -> sister -> r_cap += bottleneck;
00276 save(middle_arc);
00277 middle_arc -> r_cap -= bottleneck;
00278 for (i=middle_arc->sister->head; ; i=a->head)
00279 {
00280 a = i -> parent;
00281 if (a == TERMINAL) break;
00282 save(a);
00283 a -> r_cap += bottleneck;
00284 save(a -> sister);
00285 a -> sister -> r_cap -= bottleneck;
00286 if (!a->sister->r_cap)
00287 {
00288 set_orphan_front(i);
00289 }
00290 }
00291 save(i);
00292 i -> tr_cap -= bottleneck;
00293 if (!i->tr_cap)
00294 {
00295 set_orphan_front(i);
00296 }
00297
00298 for (i=middle_arc->head; ; i=a->head)
00299 {
00300 a = i -> parent;
00301 if (a == TERMINAL) break;
00302 save(a->sister);
00303 a -> sister -> r_cap += bottleneck;
00304 save(a);
00305 a -> r_cap -= bottleneck;
00306 if (!a->r_cap)
00307 {
00308 set_orphan_front(i);
00309 }
00310 }
00311 save(i);
00312 i -> tr_cap += bottleneck;
00313 if (!i->tr_cap)
00314 {
00315 set_orphan_front(i);
00316 }
00317
00318 flow += bottleneck;
00319 }
00320
00321
00322
00323 template <typename captype, typename tcaptype, typename flowtype>
00324 void Graph<captype,tcaptype,flowtype>::process_source_orphan(node *i)
00325 {
00326 node *j;
00327 arc *a0, *a0_min = NULL, *a;
00328 int d, d_min = INFINITE_D;
00329
00330
00331 for (a0=i->first; a0; a0=a0->next)
00332 if (a0->sister->r_cap)
00333 {
00334 j = a0 -> head;
00335 if (!j->is_sink && (a=j->parent))
00336 {
00337
00338 d = 0;
00339 while ( 1 )
00340 {
00341 if (j->TS == TIME)
00342 {
00343 d += j -> DIST;
00344 break;
00345 }
00346 a = j -> parent;
00347 d ++;
00348 if (a==TERMINAL)
00349 {
00350 j -> TS = TIME;
00351 j -> DIST = 1;
00352 break;
00353 }
00354 if (a==ORPHAN) { d = INFINITE_D; break; }
00355 j = a -> head;
00356 }
00357 if (d<INFINITE_D)
00358 {
00359 if (d<d_min)
00360 {
00361 a0_min = a0;
00362 d_min = d;
00363 }
00364
00365 for (j=a0->head; j->TS!=TIME; j=j->parent->head)
00366 {
00367 j -> TS = TIME;
00368 j -> DIST = d --;
00369 }
00370 }
00371 }
00372 }
00373
00374 if (i->parent = a0_min)
00375 {
00376 i -> TS = TIME;
00377 i -> DIST = d_min + 1;
00378 }
00379 else
00380 {
00381
00382 add_to_changed_list(i);
00383
00384
00385 for (a0=i->first; a0; a0=a0->next)
00386 {
00387 j = a0 -> head;
00388 if (!j->is_sink && (a=j->parent))
00389 {
00390 if (a0->sister->r_cap) set_active(j);
00391 if (a!=TERMINAL && a!=ORPHAN && a->head==i)
00392 {
00393 set_orphan_rear(j);
00394 }
00395 }
00396 }
00397 }
00398 }
00399
00400 template <typename captype, typename tcaptype, typename flowtype>
00401 void Graph<captype,tcaptype,flowtype>::process_sink_orphan(node *i)
00402 {
00403 node *j;
00404 arc *a0, *a0_min = NULL, *a;
00405 int d, d_min = INFINITE_D;
00406
00407
00408 for (a0=i->first; a0; a0=a0->next)
00409 if (a0->r_cap)
00410 {
00411 j = a0 -> head;
00412 if (j->is_sink && (a=j->parent))
00413 {
00414
00415 d = 0;
00416 while ( 1 )
00417 {
00418 if (j->TS == TIME)
00419 {
00420 d += j -> DIST;
00421 break;
00422 }
00423 a = j -> parent;
00424 d ++;
00425 if (a==TERMINAL)
00426 {
00427 j -> TS = TIME;
00428 j -> DIST = 1;
00429 break;
00430 }
00431 if (a==ORPHAN) { d = INFINITE_D; break; }
00432 j = a -> head;
00433 }
00434 if (d<INFINITE_D)
00435 {
00436 if (d<d_min)
00437 {
00438 a0_min = a0;
00439 d_min = d;
00440 }
00441
00442 for (j=a0->head; j->TS!=TIME; j=j->parent->head)
00443 {
00444 j -> TS = TIME;
00445 j -> DIST = d --;
00446 }
00447 }
00448 }
00449 }
00450
00451 if (i->parent = a0_min)
00452 {
00453 i -> TS = TIME;
00454 i -> DIST = d_min + 1;
00455 }
00456 else
00457 {
00458
00459 add_to_changed_list(i);
00460
00461
00462 for (a0=i->first; a0; a0=a0->next)
00463 {
00464 j = a0 -> head;
00465 if (j->is_sink && (a=j->parent))
00466 {
00467 if (a0->r_cap) set_active(j);
00468 if (a!=TERMINAL && a!=ORPHAN && a->head==i)
00469 {
00470 set_orphan_rear(j);
00471 }
00472 }
00473 }
00474 }
00475 }
00476
00477
00478
00479 template <typename captype, typename tcaptype, typename flowtype>
00480 flowtype Graph<captype,tcaptype,flowtype>::maxflow(bool reuse_trees, Block<node_id>* _changed_list)
00481 {
00482 node *i, *j, *current_node = NULL;
00483 arc *a;
00484 nodeptr *np, *np_next;
00485
00486 if (!nodeptr_block)
00487 {
00488 nodeptr_block = new DBlock<nodeptr>(NODEPTR_BLOCK_SIZE, error_function);
00489 }
00490
00491 changed_list = _changed_list;
00492 if (maxflow_iteration == 0 && reuse_trees) { if (error_function) (*error_function)("reuse_trees cannot be used in the first call to maxflow()!"); exit(1); }
00493 if (changed_list && !reuse_trees) { if (error_function) (*error_function)("changed_list cannot be used without reuse_trees!"); exit(1); }
00494
00495 if (reuse_trees) maxflow_reuse_trees_init();
00496 else maxflow_init();
00497
00498
00499 while ( 1 )
00500 {
00501
00502
00503 if ((i=current_node))
00504 {
00505 i -> next = NULL;
00506 if (!i->parent) i = NULL;
00507 }
00508 if (!i)
00509 {
00510 if (!(i = next_active())) break;
00511 }
00512
00513
00514 if (!i->is_sink)
00515 {
00516
00517 for (a=i->first; a; a=a->next)
00518 if (a->r_cap)
00519 {
00520 j = a -> head;
00521 if (!j->parent)
00522 {
00523 j -> is_sink = 0;
00524 j -> parent = a -> sister;
00525 j -> TS = i -> TS;
00526 j -> DIST = i -> DIST + 1;
00527 set_active(j);
00528 add_to_changed_list(j);
00529 }
00530 else if (j->is_sink) break;
00531 else if (j->TS <= i->TS &&
00532 j->DIST > i->DIST)
00533 {
00534
00535 j -> parent = a -> sister;
00536 j -> TS = i -> TS;
00537 j -> DIST = i -> DIST + 1;
00538 }
00539 }
00540 }
00541 else
00542 {
00543
00544 for (a=i->first; a; a=a->next)
00545 if (a->sister->r_cap)
00546 {
00547 j = a -> head;
00548 if (!j->parent)
00549 {
00550 j -> is_sink = 1;
00551 j -> parent = a -> sister;
00552 j -> TS = i -> TS;
00553 j -> DIST = i -> DIST + 1;
00554 set_active(j);
00555 add_to_changed_list(j);
00556 }
00557 else if (!j->is_sink) { a = a -> sister; break; }
00558 else if (j->TS <= i->TS &&
00559 j->DIST > i->DIST)
00560 {
00561
00562 j -> parent = a -> sister;
00563 j -> TS = i -> TS;
00564 j -> DIST = i -> DIST + 1;
00565 }
00566 }
00567 }
00568
00569 TIME ++;
00570
00571 if (a)
00572 {
00573 i -> next = i;
00574 current_node = i;
00575
00576
00577 augment(a);
00578
00579
00580
00581 while ((np=orphan_first))
00582 {
00583 np_next = np -> next;
00584 np -> next = NULL;
00585
00586 while ((np=orphan_first))
00587 {
00588 orphan_first = np -> next;
00589 i = np -> ptr;
00590 nodeptr_block -> Delete(np);
00591 if (!orphan_first) orphan_last = NULL;
00592 if (i->is_sink) process_sink_orphan(i);
00593 else process_source_orphan(i);
00594 }
00595
00596 orphan_first = np_next;
00597 }
00598
00599 }
00600 else current_node = NULL;
00601 }
00602
00603
00604 if (!reuse_trees || (maxflow_iteration % 64) == 0)
00605 {
00606 delete nodeptr_block;
00607 nodeptr_block = NULL;
00608 }
00609
00610 maxflow_iteration ++;
00611 return flow;
00612 }
00613
00614
00615
00616
00617 template <typename captype, typename tcaptype, typename flowtype>
00618 void Graph<captype,tcaptype,flowtype>::test_consistency(node* current_node)
00619 {
00620 node *i;
00621 arc *a;
00622 int r;
00623 int num1 = 0, num2 = 0;
00624
00625
00626 for (i=nodes; i<node_last; i++)
00627 {
00628 if (i->next || i==current_node) num1 ++;
00629 }
00630 for (r=0; r<3; r++)
00631 {
00632 i = (r == 2) ? current_node : queue_first[r];
00633 if (i)
00634 for ( ; ; i=i->next)
00635 {
00636 num2 ++;
00637 if (i->next == i)
00638 {
00639 if (r<2){ assert(i == queue_last[r]);
00640 }else assert(i == current_node);
00641 break;
00642 }
00643 }
00644 }
00645 assert(num1 == num2);
00646
00647 for (i=nodes; i<node_last; i++)
00648 {
00649
00650 if (i->parent == NULL) {}
00651 else if (i->parent == ORPHAN) {}
00652 else if (i->parent == TERMINAL)
00653 {
00654 if (!i->is_sink){ assert(i->tr_cap > 0);
00655 }else assert(i->tr_cap < 0);
00656 }
00657 else
00658 {
00659 if (!i->is_sink){ assert (i->parent->sister->r_cap > 0);
00660 }else assert (i->parent->r_cap > 0);
00661 }
00662
00663
00664 if (i->parent && !i->next)
00665 {
00666 if (!i->is_sink)
00667 {
00668 assert(i->tr_cap >= 0);
00669 for (a=i->first; a; a=a->next)
00670 {
00671 if (a->r_cap > 0) assert(a->head->parent && !a->head->is_sink);
00672 }
00673 }
00674 else
00675 {
00676 assert(i->tr_cap <= 0);
00677 for (a=i->first; a; a=a->next)
00678 {
00679 if (a->sister->r_cap > 0) assert(a->head->parent && a->head->is_sink);
00680 }
00681 }
00682 }
00683
00684 if (i->parent && i->parent!=ORPHAN && i->parent!=TERMINAL)
00685 {
00686 assert(i->TS <= i->parent->head->TS);
00687 if (i->TS == i->parent->head->TS) assert(i->DIST > i->parent->head->DIST);
00688 }
00689 }
00690 }